Probing the energy bands of a Bose-Einstein condensate in an optical lattice 
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We simulate three experimental methods which could be realized in the laboratory to probe 
s ! ■ the band excitation energies and the momentum distribution of a Bose-Einstein condensate inside 

' an optical lattice. The values of the excitation energies obtained in these different methods agree 

within the accuracy of the simulation. The meaning of the results in terms of density and phase 
■ deformations is tested by studying the relaxation of a phase-modulated condensate towards the 

ground state. 
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(N ! PACS numbers: 67.40.Db, 03.75.-b, 32.80.Pj 

A periodic potential can be imposed on a cold atomic gas from the shift of the atomic ground state under illumination 
by a detuned laser standing wave Q. Experiments on ultracold atoms inside such an optical lattice have revealed 
dynamical behaviours which are well known from band theory for electrons in solids, i.e. Wannier-Stark ladders [^), 
Bloch oscillations || and Landau-Zener tunneling [Q . Especially exciting perspectives are offered by confinement of an 
+i> \ atomic Bose-Einstein condensate (BEC) inside an optical lattice, in regard to the generation of coherent matter-wave 
pulses laser cooling ||, quantum computing Jt| and more generally BEC quantum-state engineering. 

Theoretical and numerical studies of a BEC in a periodic potential have been based on the Bogolubov-de Gennes 
equations ||- |llj and on the Gross-Pitaevskii equation (CPE) [^2)- These studies have shown that, assuming 
phase coherence in the ground state, a BEC subject to a periodic optical potential and to a constant external force 
performs Bloch oscillations as if it were a quasi-particle inside the lowest energy band. This is consistent with the 
observed emission of coherent matter pulses from a vertical array of optical traps under the force of gravity Q . 

In this Letter we report the first detailed theoretical analysis of the observability of the BEC energy bands under 
currently attainable experimental conditions. To this end we use the GPE to simulate three different pump-probe 
experimental methods, which have already been realized in the laboratory to measure Bloch oscillations of ultracold 
Bose atoms || and to study the shape deformation modes of a BEC in harmonic traps |l^,^8|. In interpreting our 
results it will be useful to recall that (i) the lowest band arises from pure phase modulations of the BEC with given 
q- vector in the Brillouin zone, and (ii) the higher bands are associated with BEC density profiles having the symmetry 
of the isolated- well states, accompanied by modulations of the phase (see especially Ref. Jll|], where these properties 
are derived in the Wannier representation for the wave functions). 
C^l ■ We consider a dilute BEC in the external potential U®[1 — exp(— r 2 /rf b ) cos 2 (27rz/A)], with U? being the well depth, 
A and rn, the laser-beam wavelength and waist, and d = A/2 the lattice period. Its dynamics is accurately described 
by a one-dimensional (ID) model, if the mean-field interactions are renormalized to reproduce the correct 3D value 
of the chemical potential fl4|| . We thus adopt the ID GPE for the BEC wave function ^>(z,t), 
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In Eq. (Q) U (z, t) includes the ID lattice potential Ui(z) and the mean-field interactions, i.e. U(z, t) = sin 2 (Trz/d) + 
4Trh 2 jap\ 1 i'(z, t)\ 2 /M, with a the scattering length, 7 the renormalization factor, p the number of particles per lattice 
well and M the atom mass. We study three different forms for the pump potential U p (z, t), as will be specified below. 
We adopt the system parameters from the experiment on 87 Rb [^) (a — 110 Bohr radii, A = 850 nm and — 1.4 Er 
with Er — h 2 /8Md 2 ) and take p = 2500 atoms per well. 

On the technical side, we handle the heavy simulations needed to obtain significant results on fine spectral struc- 
tures by a well tested explicit-time-marching algorithm (T^|l{|. We first determine the ground state by numerically 
propagating Eq. ([j]) in imaginary time and then insert it as the initial condition for evolution in real time. By the 
same algorithm we preliminarly calculate the band energies E n (q) from the relaxation of a statically deformed BEC 
(see below), n being the band index and q the reduced wave vector in the Brillouin zone. The size of the simulation 
box is 700 wells, as needed for a three-digits accuracy in the ground-state energy (p = 0.695 En). We guard against 
unwanted localized excitations at the boundaries by letting the density profile vanish over a length scale which is much 
larger than the healing length £ = (Snap)^ 1 ^ 2 ~ 1.4 d. Typically, a grid contains 5 • 10 4 points and the simulation is 
carried out up to final times ranging from 13 to 32 ms. In the dynamical simulations the data are stored every 20 ps 
and used to obtain the real and imaginary parts of the Fourier transform of the wave function (^(z, to), say). 



1 



We turn to present our numerical results. In the first method, which we familiarly refer to as kicking, we impart a 
velocity v — Tiq/M to the BEC at time t = by imposing a phase hqz/M on the ground-state wave function. The 
dynamics of the BEC is then monitored with U p = 0. We show in Fig. 1 the spectra uS) at any position z in 

the bulk of the BEC, as functions of {Two — [x)/Er for four values of q. At each q the strong peak on the left refers 
to the lowest band (n = 0, say) and the other two peaks belong to the bands n — 1 and n — 2. The eigenfrequencies 
revealed by Fig. 1 are the same as those appearing in 3^(2, u/). The density spectra \^(z, w)| 2 show instead a strong 
peak centered at energy /i for all values of q and two further peaks corresponding to the n — 1 and n — 2 peaks in 
Fig. 1, the n — peak of Fig. 1 being absent. 
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FIG. 1. Kicking method: bulk spectrum $t&(z,u) as a function of (hui — fij/En for four values of q in the Brillouin zone. 
From bottom to top: q = 0.25, 0.50, 0.75 and 1.0 n/d. 

Thus, a measurement of the density |\&(Az, t)\ 2 of a portion Az of the BEC by a technique of imaging after ballistic 
expansion [ fL7| could reveal excitations in higher bands, but should have no access to those in the lowest band. This is 
in accord with the aforementioned theoretical result ]Tl[ | showing that the excitations in the n — band correspond to 
pure phase modulations. Instead, the momentum distribution n(p,t) = \J dz^(z, t) exp(ipz)\ 2 contains information 
on the spatial variations of the BEC phase through the average velocity. The spectra shown in Fig. 1 could therefore be 
experimentally accessible through measurements of the momentum distribution n(p, t) as a function of the observation 
time after imparting an initial velocity (see also the discussion below). 

Turning to the second simulational method, that we refer to as shaking, we modulate the BEC in space and time 
by parametrically driving the lattice potential. To this end we set U p (z, t) = alii (z) cos(qz — fit) for t < and U p = 
otherwise. We choose a small pump amplitude (a = 0.15), match the pump frequency Q near resonance and tune the 
drive-time t& over several excitation periods, the time unit being Tr = h/ Er — 0.32 ms. The density is then recorded 
at times t^>td- 

Fig. 2 displays the bulk density spectrum ^(zjw)! 2 for the same q-values as in Fig. 1. The spectrum contains 
a peak centered at energy /x for all values of q and a strong peak corresponding to the nearly resonant drive, at a 
frequency changing with q in accord with the n = peaks in Fig. 1. The other structures in Fig. 2 can be identified 
as harmonics of the fundamental excitation energy and combinations of them. 

Thus a near-resonance measurement of the density \^f(z, t)\ 2 in the shaking method may reveal the excitations in 
all bands and in particular those in the n = band. Shaking at q = coherently drives the bound states in each well 
and creates pure density fluctuations. Phase modulations are also triggered at q ^ 0, due to the spatial periodicity of 
the BEC: this is at variance from the shape-deformation modes in harmonic traps ]l8|| . The energies of these phase 
and density excitations are located at the poles of the response function to the external drive and therefore show up 
as peaks in \9(z, ^)| 2 - 
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FIG. 2. Shaking method: bulk density |*(,z,tj)| 2 as a function of (hui — h)/Er for the same g-values as in Fig. 1. 

Before quantitatively comparing the results obtained by the above dynamical methods, we pause to present the 
band structure that we have obtained by a static deformation method. In this method we imprint a phase modulation 
on the BEC ground state so as to generate a state having some overlap with a Bloch state of quasi-momentum q. We 
then propagate this modulated BEC in imaginary time and monitor the average energy < E > as the BEC returns 
to its ground state. Fig. 3 shows < E > vs. UTi/Er for the same q- values as in Fig. 1. Two plateaus are met during 
this evolution: the first lies at energy E (q) and the second at energy \x. The time scales for the appearance of the 
plateaus correspond to interband and intraband relaxation, respectively. Similar results are obtained for the higher 
bands E n (q) from phase modulations having wave vector q + n'ir/d outside the first Brillouin zone (with n' — n for 
even n and n' = n + 1 for odd n) . 
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FIG. 3. Static deformation method: average energy < E > vs. imaginary time ith/En for a modulated BEC with the same 
q values as in Figs. 1 and 2. Higher plateaus correpond to higher values of q in the lowest band. 



3 



Fig. 4 reports the energy bands E n (q) that we have obtained by the three methods presented above, for n = 0, 1 
and 2. The three methods yield the same excitation energies within the simulation accuracy which is reflected by the 
size of the symbols (squares for kicking, circles for shaking and triangles for the static method) in the main body of 
the figure. A gap has opened at q = ir/d between the n = and n = 1 bands, while there is no gap at q = between 
the n = 1 and n — 2 bands: this is a well known result for ID lattice potentials having a doubled period p0| . We 
also remark that the results shown in Fig. 4 for the n = 1 and n — 2 bands do not differ significantly from those that 
we obtain for a non-interacting BEC: indeed, for the present system parameters we have Uf ~ 2.12/i, implying that 
mean-field interactions are weak on the energy scale of the lattice. 
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FIG. 4. The first three energy bands in the Brillouin zone < q < n/d as obtained by the methods reported in Figs. 1-3 
(squares: kicking; circles: shaking; triangles: static method). The inset zooms on the lowest band. The size of the error bars 
reflects a two-digit accuracy in energy differences. Solid line: quadratic dispersion with the bare atom mass. 

The inset in Fig. 4 shows an enlarged view of the excitation energies in the lowest band up to qd/n = 0.5, where 
the band starts bending over. On this scale the estimated error bars of the simulation become visible and well within 
them the results for the present interacting BEC agree with those for the non-interacting one, except perhaps at 
qd/n = 0.1. The calculated band energies are compatible with the free-particle dispersion relation h 2 q 2 /2M as shown 
by the solid line in the inset. High resolution and fine tuning of the strength of the self- interactions relative to the 
height of the lattice confinement would clearly be needed to experimentally reveal the phonon-like dispersion relation 
at long wavelengths which is predicted by gapless theories p],[y] ■ From the Bogolubov dispersion relation this linear 
behaviour may be expected to become visible below q w v2/£. 

Finally, we return to discuss the momentum distribution n(p, t) for a BEC in an optical lattice. We simulate a 
feasible experimental method for its measurement by driving the BEC with a force F for a variable amount of time 
td. The force is taken as a positive constant in the range < td < T^/2 and as a negative constant in the range 
Tg/2 < id < Tg, with Tg = h/\F\d being the period of Bloch oscillations. We set U p (z,t) — —Fz for t < td and 
Up = otherwise, and choose the value of F to correspond to a particle acceleration of 85 cm/ s 2 as in the experiments 
of Ben Dahan et al. ||. The Brillouin zone is explored in this method according to qd/ir — 2sgn(F)td/TB- 

Fig. 5 shows n(p, t) at times < td < Tg covering one full Bloch oscillation. The peak seen at td = drifts to 
the right up to the first zone boundary at td = Tb/2, while a second peak emerges in the second zone and moves to 
the left entering the first zone at later times. The sharpness of the peaks shows that the BEC is behaving as if it 
were a quasi-particle reflected back and forth at the Bragg planes. The inset in Fig. 5 shows that the BEC average 
velocity < v >, as obtained from the momentum distribution after backfolding of the second zone, is consistent with 
the semiclassical average velocity as calculated from dEo(q)/dq using the results in Fig. 4. 
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FIG. 5. Force-driving method: momentum distribution (in arbitrary units) at various times id, equally spaced by Ts/16 up 
to Tb/2 and then by Tb/8 up to Tg. Inset: points give the BEC average velocity < v > (in units of vr — hn/Md) as determined 
from the momentum distribution after backfolding of the second Brillouin zone, while the line shows the semiclassical average 
velocity. 

In conclusion, we have simulated three different experimental methods to probe the energy bands of a condensate 
in an optical lattice and tested their meaning in terms of the density and phase deformations associated with excited 
band states. The kicking and the shaking methods do not require any special modelling of the observables in an actual 
experiment. The force-driving method is especially suited for a measurement of the BEC momentum distribution 
in an optical lattice. A study of the observability of the expected linear dispersion at small g-values would require 
different system parameters than the present ones and in particular stronger mean-field interactions and larger lattice 
constants. Theoretical and numerical work is in progress to explore these possibilities. 

One of us (MLC) thanks Dr J.H. Miillcr for enlightening discussions on experimental issues. Part of this work was 
done using the computational facilities of Cineca. We acknowledge support from MURST through PRIN2000. 
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